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The biological membrane, which compartmentalizes the cell and 
exhibit wide variety of macroscopic shapes of varying morphology 
and topology.s A systematic understanding of the relation of mem- 
brane shapes to composition, external field, environmental condi- 
tions etc. have important biological relevance. Here we review 
the triangulated surface model, used in the macroscopic simula- 
tion of membranes and the associated Monte Carlo (DTMC) meth- 
ods. New techniques to calculate surface quantifiers, that will facil- 
itate the study of additional in-plane orientational degrees of free- 
dom, has been introduced. The mere presence of a polar and ne- 
matic fields in the ordered phase drives the ground state confor- 
mations of the membrane to a cylinder and tetrahedron respec- 
tively. 



its organelles, 




Introduction 

The biomembrane, that define the shapes of the cells and organelles that they enclose, 
are a few nanometer thick, with their lateral size extending from ten's of nanometers up 
to microns ^ . These membranes can organize itself into a wide array of shapes, from 
relatively smooth shapes, as in plasma membranes, to complex cisternae, as in the ER 
and Golgi. Understanding the properties of membranes and their interaction with other 
macromolecules, at multiple scales, is absolutely essential in modeling most biological 
processes. However, the complexity arising from the chemical diversity of membrane 
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components hinders experimental investigations and modeling at the molecular scale. In 

El E3 

the classical model proposed by Canham and Helfrich , the homogeneous membrane, 
with its components having a diffusive degree of freedom, represents a two dimensional 
fluid surface with conformational energy, 

■W sur = J dS ^H 2 + a + n G R^ + J ApdV. (1) 

The phenomenological parameters, k, a, k g and Ap, respectively the bending rigidity, 
surface tension, Gaussian rigidity and osmotic pressure difference, depend on the chemi- 
cal composition of the membrane. H and R are the mean and Gaussian curvature of the 
surface. For an excellent review on the applicability of this and other related models see 
the article by Seifert. 

The anisotropy in bending modulus, arising from the presence of embedding proteins 
and nanostructured membrane domains, have been found to influence membrane shapes. 
For instance, the presence of surface membrane proteins containing BAR domains, like 
Nexin, Dynamin, Caveolin etc., have been shown to induce a wide array of shapes to the 
membrane. Coatomers, Dynamin and F-BAR domain containing proteins stabilize 
tubular membrane structures, while proteins like Caveolin, I-BAR and some class of tox- 
ins like the Shiga toxin form caveola, that are tubular structures drawn into cytoplasm. 
Protein induced membrane curvature is believed to be one of the possible mechanisms con- 
trolling membrane morphology. This idea is further supported by molecular simulations 

where proteins oriented in an ordered phase, have been shown to induce a spontaneous 
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curvature to the membrane. The protein distribution and their cooperative effect 
on membrane morphology has been theoretically investigated by extending the existing 
membrane models. However, the realm of these theoretical models are restricted to 
axisymmetric membrane shapes and small deviations around it. 



Triangulated surface model for membrane 

In this approach, the vesicle, a two dimensional closed surface, embedded in three dimen- 
sional space, is discretized into a set of, interconnected, self avoiding, triangular plaque- 
ttes. Each of these plaquette correspond to a flat bilayer patch of the membrane. The T 
plaquettes, constituting L independent links, intersect at iV vertex points. The topology 
of the closed discretized surface, is defined by the Euler characteristics x — N + T — L. 

The implicit discretization of the elastic energy of the membrane in EquationQ is 
given by 

* = Afe E E i 1 - ^[/i( e )]-^[/2(e)]} + Ap V. (2) 

« {e}v 



As illustrated in Figure 1(a) , {e(v)} denote the set of links, a vertex v makes in its one 



ring neighbourhood. N[fi(e)] and N[f 2 (e)} are respectively the normals to faces fi and 
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Figure 1: (a) A surface patch showing a one ring neighbourhood around vertex v. N(v) 
represents the normal to the tangent plane at v. R{e) represents the vector along an edge 
e, while N(e) is its normal, (b) the signed dihedral angle $(e) between faces, fi and /2, 
sharing an edge e. 



f2 sharing an edge e. V is the volume enclosed by the surface. The bending rigidities in 
Equations^ and ([I]) are related to each other as, = \/3k for a sphere and At = 2/t/v / 3 
for a cylinder. EU The triangulated mesh can be evolved using various techniques suiting 
the nature of study. We will, in this article, concentrate on the equilibrium properties 
of the membrane using Monte Carlo methods for mesh evolution. Details pertaining to 
simulations can be found in literature, also see chapter by Gompper and Kroll^ 
for an excellent review on the application of triangulated surface techniques to various 
classes of problems. 



In-plane orientational order 

The additional degrees of freedom, that arise from the anisotropy induced by the proteins, 
peptides, cytoskeletons and from tilt of lipids in the bilayer, can be quantified by a p- 
atic vector living on the membrane surface. For example, the nematic order is defined 
as a 2— atic vector field, with a tt rotational symmetry, in the local Darboux frame, 
constructed from the principal directions and the vertex normal at each vertex of the 
mesh. This means we need to calculate the principal curvatures and directions, unlike 
in EquationQ where the squared mean curvature is directly approximated from the 
face normals. We will outline, in this article, the basic steps to compute these surface 
geometrical quantifiers. ^ 
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Computing surface quantifiers 

From the neighborhood of a vertex v, the normal to the surface at the vertex is calculated 



as . 



N(v) 



{fh 



I 



{fh 



where A(f) is the area and N(f) is the unit normal at face /. The weight factor Q is 
chosen to be proportional to the area A(f). The normals of the faces sharing an edge e, 
can be used to approximate the normal to the edge as, 



N(e)= fi\Ji{e)] + ft[f 2 (e)) / N[h(e)} + N[f 2 (e)} 



To quantify the curvature, that arises from the faces being non-planar, we define 



H(e) 



R(e) 



cos 
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along the direction B(e) 



N(e) x R(e). mm $(e 



signed dihedral angle between the faces, /i(e) and /2(e), 



is the 

taking a value of it when 



the faces are coplanar 



R(e) 



is the length of the edge. The discretized "shape 



operator", that represents the curvature tensor and its orientation at e is given by 



S e (e) = H{e) 



B(e)®B(e) 



. It should be noted that S e (e) is defined at every edge con- 
taining vertex v, which can collectively be represented by {S e (e)}. These individual shape 
operators determine the curvature at v. The projection operator, P(i>) = l — N(v)(&N(v), 
projects S e (e) into the tangent plane at vertex v. Details of the surface quantifiers at 
v are contained in the vertex shape operator, constructed as a weighted sum of these 
projections, given by, 



SJv) 



A(v) 



^W{e)P(v)i$Jie)P(v). 

{e}v 



(3) 



A[v) = ^(/)/3 is the average surface area around v, while the weight factor for an 

edge is calculated as W(e) = N(v) ■ N(e). S v (v) is constructed in the global cartesian 

frame and a Householder transformation, rotates the operator in Equation(j3]) to 
its tangent frame and the resulting matrix is a 2 x 2 minor. This transformation gives 
a computationally efficient route to calculate the eigenspectrum of Equation ([3]), when 
compared to using standard numerical techniques. The eigenvalues and eigendirections 
of S v (f) are respectively the principal curvatures(ci & c 2 ) and directions(ei & e%). The 
scalar invariants are the mean curvature, H — (ci + oi)j1 and the Gaussian curvature, 
R = C\C2- We can now define the Darboux frame with basis vectors [ei(i>), e 2 iy ), N(v)] 
and the in-plane nematic orientation, in this frame, to be fh = cos(p(v)ei + smip{v)e2- 
(f(v) is the angle subtended by the nematic orientation with respect to maximum principal 
direction e\. 
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Parallel transport on a discrete mesh 



In order to compare the orientation of two distant in-plane vectors on the surface, it is 
necessary to perform a parallel transport of the vectors on the discretized surface. In 
practice, we need only to define the parallel transport between neighboring vertices, i.e. 
a transformation rh(v ) — > T_(v, v )rh(v), which brings m(v) correctly into the tangent 
plane of the vertex v , so that its angle with respect to the geodesic connecting v and 
v is preserved. If r(v,v ) is the unit vector connecting a vertex v to its neighbor v and 
((v)=P(v)r(v, v) and ((v')= P_(v')r(v\v) are its projection on to the tangent planes at v 
and v ; then our best estimate for the directions of the geodesic connecting them, are the 
unit vectors C(v), ((v ). The decomposition of m(v ) along the orientation of the geodesic 
and its perpendicular in the tangent plane of v is thus: 



m(v 



) = [m(v) ■ C(v)} ((v) + [m(v) ■ [N(v) x ((v)}} [n(v) x ((v)\ . (4) 



Parallelism now demand that these coordinates, with respect to the geodesic orientation, 
are the same in the tangent plane of v , therefore: 



,«>(») = {m(«) • ((»)}((«') + (m(«) • [JV(d) x ((«)]} [n(v') x C(v')] . (5) 



T(v 



This parallel transport operation allow us to define the angle (p(v,v) between vectors in 
the tangent plane at neighboring vertices, and in turn their cosine and sine as: 



cos(</?(i>, v )) 
sm({p(v, v )) 



™( v ) ' E(i>, v )rh(v); 
N(v ) x rh(v ) ■ T(v, v )rh(v) 



(6) 



We can now define the lattice models with interaction between in-plane orientational 
field, e.g., the XY-model for a polar order {p = 1) on a random surface: 



"Hi-atic = — ^- ^2 cos[<p(y,v')] 

(vv') 



(7) 



or the Lebwohl-Lasher model for nematic order(p = 2) on a random surface: 

"H2-atic 



t22 \ ^ 

{vv ) 



^cos 2 (v?(^,t/)) - X - 



(8) 



We will, in general, represent the self interaction energy of a p-atic field by J^- a tic- 
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Figure 2: A Monte Carlo move involving the in-plane field, fa. The membrane morphol- 
ogy remains unchanged while the field orientation is changed from ip to tp . e\ and £2 are 
the principal directions. 



Frustrated in-plane order and membrane morphology 

The curvature inducing properties of surface proteins, can be captured in a phenomeno- 
logical model as an additional term to the RHS of Equation ([I]), given by, 

M> NM = $>(«) {^l H \\( v ) - 4] 2 + - 4] 2 } • (9) 

V 

K\\ and k± are the directional bending rigidities parallel (or antiparallel) and perpendicular 
to the field m. These phenomenological parameters, along with the directional sponta- 
neous curvatures cjj and qj- , quantify the interaction strength of the protein with the 
membrane. The directional curvatures of the membrane at a vertex v, calculated using 
Gauss formula, ^ are Hn(v) = C\(v) cos 2 (p(v) + c%(v) sin 2 <^(t>) and H± = 02(f) cos 2 <p(v) + 
c\{v) sin 2 <p(v). The additional Monte Carlo move involving the in-plane order is given in 

Figure 2 ^ 

The total energy of the field decorated surface, of spherical topology (genus, g — 0), is 
J^tot = J^lur + '^NM + ^p-atic- In the ground state of J#tot, the p-atic field is frustrated by 
the topology of the embedding surface and has in it XP defects, each of topological charge 
1/p. In addition to the bending rigidity, the equilibrium conformation of the surface now 
depends on the number and strength of the defects. This implicit coupling of membrane 
morphology to the in-plane order and the resulting membrane shapes has been theoreti- 
cally investigated in earlier works. The equilibrium conformations of the membrane 
in our simulations, for k\\ — k± — 0, are in excellent agreement with these predictions. 
Figure [3(a)| shows the predicted ellipsoidal membrane for a polar field (p = 1), while the 



tetrahedral shape, characteristic of a nematic field, is shown in Figure 3(b) In case of the 
polar fields the two +1 antipodal vortices are positioned at the ellipsoidal caps whereas 
the four +1/2 disclinations are at the vertices of the tetrahedron for the nematic field. 

In the case of a fluid membranes, the shape change induced by the ordered in-plane 

I221 

fields was investigated earlier by us. It has been shown that the entropy dominated 
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+ 1 vortex 
at the poles 




Ellipsoidal membrane 



(a) 



+1/2 disclinations at the 
vertices of the tetrahedron 




u \o \\A\ 



A\f 



Tetrahedral membrane 



(b) 



Figure 3: A surface of spherical topology with (a) polar order (p = 1, en = 10) : Two 
+ 1 defects sit at the pole and an otherwise spherical membrane elongates into a cylinder 
(b) nematic order (p = 2, 622 = 10) : Four +1/2 disclinations remodel the surface into a 
tetrahedron. Simulations were performed with (3 = 1/fceT = 5.0 and n = 10 

branched shapes, the equilibrium shapes of flexible surfaces(/t = 0), are cut off by or- 
der induced membrane stiffening. These observations add support to the hypothesis of 
membrane shape stabilization by protein-lipid interactions. ^ Further, non zero contri- 
butions from the directional rigidities and spontaneous curvatures drive the membrane 
into a wide array of shapes encompassing the biologically relevant tubes, discs, caveolae 
and branches. ^ In this article, using a polar in-plane field, we demonstrate the emer- 
gence of both axisymmetric and non axisymmetric membrane conformations as in Figure 
[1J The polar field interaction with the membrane is thread like (k± = 0) and hence 
deforms the membrane only along its long axis. 



The biconcave and twisted disc shapes shown in Figure 4(a) and 4(b)| are, respectively, 
the ground state configurations for Cq = 0.25 and C = 0.75. The in-plane polar fields 
in Figure [4] have similar defect structures, two +1 vortices in this case, and hence have 
comparable energies of M'x-^ic- Thus, the equilibrium membrane deformation is driven by 
the competition between the elastic energy (Equation Q) and field-membrane interaction 
energy (Equation (|9])). The minimum energy configuration of EquationQ is a minimal 
surface(i/ = 0) while that of Equation(|9| corresponds to a surface with the field oriented 
along the maximum principal direction e 1; and the corresponding curvature c\ = Cq. 
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Figure 4: Conformations of a fluid membrane with in-plane polar(jo = 1) field for (a) 
Cq = 0.25 and (b) Cq = 0.75. Values of other parameters are k — 10, /3 — 5.0, Ku = 5 
and en = 10. 

The observed trend in the values of the principal curvature, as a function of Cq, are in 
agreement with the above arguments ^ . 



Conclusion 

We have reviewed the basic techniques of dynamically triangulated surfaces used as model 
for fluid membranes with in-plane order in the macroscopic limit. The anisotropic mem- 
brane constituents are modeled as in-plane vector fields and surface rendering techniques 
suitable to study these systems have been introduced. The presence of an ordered in- 
plane p-atic field remodels the membrane into an ellipsoid(p = 1) and cylinder (p = 2). 
Introducing explicit field-membrane interactions widens the spectrum of the resulting 
shapes and can serve as an useful way of understanding the cooperative effect of protein 
lipid interactions. 
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